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ABSTRACT 


Traditional  analyses  of  optical  scintillation  usually  invoke  the  Rytov  approximation  and  the 
associated  assumptions  regarding  weak  phase  perturbations  throughout  the  propagation  channel. 
Associated  with  that  approach  is  the  prediction  that  the  statistics  of  scintillation  will  be  log-normal.  This 
report,  which  uses  the  Fresnel  propagation  function  as  a  basis,  challenges  many  of  the  standard 
descriptions  of  far-field  irradiance  fluctuations.  A  key  finding  of  this  study  is  the  importance  of  the  ratio 
between  the  transmitted  beam  diameter  and  the  turbulence  coherence  diameter.  Turbulence-induced  beam 
jitter  is  found  to  be  a  dominant  effect  when  this  radius  is  close  to  unity,  and  the  relationship  between 
pointing  error  and  scintillation  is  examined  in  detail.  As  a  result  of  this  work,  models  for  the  mean  Strehl, 
scintillation  index,  and  signal-fade  distributions  have  been  developed.  Included  in  this  document  are  a  set 
of  detailed  comparisons  between  the  new  theoretical  models  and  numerical  results  derived  from  a 
beamwave  propagation  code  developed  by  Lincoln  Laboratory. 
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1.  INTRODUCTION 


Scintillation  effects  due  to  electromagnetic  wave  propagation  through  an  inhomogeneous  medium 
are  generally  treated  as  a  weak  perturbation  of  the  vacuum  solution  to  Maxwell’s  electromagnetic  field 
equations.  When  this  approach  is  combined  with  a  statistical  description  of  turbulence-induced  refractive 
index  fluctuations,  one  can  accurately  model  the  behavior  of  far-field  irradiance  fluctuations  when  the 
path-integrated  wavefront  distortion  is  small  compared  to  the  optical  wavelength.  Although  many 
investigators  have  also  applied  this  analytical  approach  to  strong  turbulence  conditions,  the  results  of  this 
study  demonstrate  that  better  descriptions  of  the  scintillation  phenomenon  result  when  the  propagation 
expressions  are  modified  to  permit  a  wider  range  of  validity. 

The  need  for  a  nontraditional  description  of  optical  scintillation  is  immediately  apparent  to  anyone 
who  has  carefully  studied  beam  images  generated  by  numerical  simulations  or  experiments  in  which  laser 
beams  are  propagated  through  an  atmospheric  channel.  The  compositional  evolution  of  the  short-exposure 
beam  profile  is  dramatic  as  the  strength  of  the  field  perturbations  increases,  and  therefore  no  single 
description  adequately  characterizes  the  entire  range  of  behavior.  This  report  identifies  three  distinct 
performance  regimes  and  derives  separate  statistical  models  for  each. 

This  study  is  primarily  focused  on  the  development  of  better  link-budget  models  for  free-space 
optical  communication  systems,  but  the  derived  results  are  applicable  to  a  wide  range  of  propagation 
scenarios.  The  statistical  descriptions  developed  thus  far  only  relate  to  uncorrelated  path  realizations. 
Time-dependent  behavior  will  be  addressed  in  this  next  phase  of  this  investigation. 
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2.  THE  RYTOV  APPROXIMATION  FOR  WEAK  PERTURBATIONS 


To  establish  an  historical  perspective  for  the  new  work  discussed  in  subsequent  sections,  it  will  be 
useful  to  briefly  review  the  standard  assumptions  and  predictions  that  derive  from  the  application  of  the 
Rytov  perturbation  theory.  The  Rytov  approximation  is  described  in  a  variety  of  texts;  the  following 
derivation  is  extracted  from  Goodman  [1]. 

At  any  point  in  space  the  optical  field,  U,  must  satisfy  Maxwell’s  equation 

V2U  +  ^-U  =  0  ,  (2-1) 

c 

where  co  is  the  angular  frequency,  n  is  the  refractive  index,  and  c  is  the  speed  of  light.  For  an 
inhomogeneous  medium,  the  refractive  index  can  be  divided  into  a  constant  term,  no ,  and  a  perturbation, 
/7 / ,  that  will  vary  in  space  and  time.  The  expression  for  the  perturbed  field,  U  ,  is 

V2{U  +  U£)  +  ^-(n+n£f  (U  +  Us)-0  ,  (2-2) 

c 

where  U  and  Ue  are  the  unperturbed  and  perturbation  field  components,  respectively.  For  atmospheric 
propagation,  n  can  be  set  to  unity  without  loss  of  generality.  If  n£  is  small  compared  to  unity,  then  the 
following  expression  will  be  valid 

V2  U£  +  k2  U e  =  —2  k2  n£  U  ,  (2-3) 

where  k  =  coj c  =  2/r/ A  and  A  is  the  optical  wavelength.  For  plane-wave  propagation  through  frozen 
turbulence,  the  solution  to  Eq.  (2-3)  can  be  expressed  as  a  three-dimensional  coherent  integral  over  the 
entire  scattering  volume 


1 

"  exp] 

[  i  k\ 

\r  - 

-4) 

4  K 

J 

|  r  - 

— f  i 

r 

v 


2  k2ns(f')U(i 


/  3 

a  r 


(2-4) 


Goodman  notes  that  the  effective  scattering  angle  due  to  refractive-index  fluctuations  is  very  small,  in 
which  case  the  Fresnel  approximation  can  be  applied.  The  separation  of  the  axial  and  radial  components 
in  the  free-space  Green’s  function  yields  the  following  form  for  beam  propagation  along  the  positive  z 
axis 
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The  relationship  given  in  Eq.  (2-2)  is  referred  to  as  the  Born  approximation,  which  assumes  that  the 
net  field  value  can  be  represented  as  the  linear  sum  of  the  unperturbed  field  and  a  set  of  weak  perturbation 
components.  An  alternative  approach,  known  as  Rytov  approximation,  defines  the  complex  phase, 
y/  -  x  +  ify ,  such  that 

0  =  exp [y/]  =  exp {ys  +  y/€)  .  (2-6) 

In  this  approximation  the  phase  perturbation  components  are  assumed  to  be  additive.  For  weak 
fluctuations 


Ve  =  ln 
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so  that  the  volume  integral  can  be  written  as 


(2-7) 
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The  real  part  of  this  expression  specifies  the  log-amplitude  component  of  the  field  variations  due  to 
inhomogeneities  in  the  refractive-index.  The  Rytov  approximation  is  more  frequently  used  than  the  Born 
approximation  and  is  thought  to  be  more  accurate  for  strong  turbulence  conditions. 

Statistical  investigations  of  scintillation  due  to  refractive  index  fluctuations  in  the  atmosphere 
usually  begin  with  an  expression  similar  to  that  given  in  Eq.  (2-8)  and  proceed  with  an  analysis  based  on 
the  spatial  frequency  properties  of  the  perturbation  parameter,  Analytical  investigations  of  this  type 
date  back  to  the  work  of  Kolmogorov  [2]  and  Tatarski  [3],  and  have  since  been  continued  by  a  number  of 
investigators  including  Ishimaru  [4],  Strohbehn  [5],  Fried  [6],  and  Yura  [7].  The  results  of  these 
investigations  are  typically  expressed  in  terms  of  complex  multidimensional  integrals,  which  are  usually 
evaluated  numerically.  A  more  recent  technique  based  on  Mellin  transforms  has  been  developed  by 
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Sasiela  [8],  which  allows  the  solutions  to  many  important  problems  to  be  represented  in  closed  form.  To 
simplify  the  mathematics,  most  analytical  calculations  are  limited  to  two  propagation  geometries — 
irradiation  from  a  point  source  or  an  infinite  plane  wave. 

The  computational  approaches  just  described  have  been  applied  successfully  to  a  wide  variety  of 
practical  problems,  including  the  design  of  adaptive-optics  correction  systems  for  astronomical  imaging 
and  link  budget  estimation  for  free-space  optical  communications.  The  primary  characterization  of  the 
atmospheric  channel  is  the  Cn  profile,  which  quantifies  the  strength  of  the  refractive  index  fluctuations. 
For  many  system  constructs  the  mean  receiver  irradiance  and  the  irradiance  variance  can  be  related  to  the 
low-order  moments  of  the  Cn  profile,  i.e., 

L 

Mm  =  JV(Z)Z"'  CHz)dz  ’  <2'9) 

0 

where  L  is  the  distance  between  the  transmitter  and  receiver,  and  /describes  the  beam  dimensions.  The 
standard  result  for  the  variance  of  x  for  point-source  radiation  propagating  from  ground  to  a  receiver  in 
space  at  range  L  is 

~  0-563  &7//6  /v5/6  ,  (2-10) 

and  for  an  infinite  plane  wave  traversing  the  same  path 

a2x  »  0.563  k1/b  /r0  LS/6  .  (2-11) 

Based  on  the  central  limit  theorem  it  is  usually  argued  that  both  ^and  ^are  Gaussian  random  variables,  in 
which  case  Eq.  (2-7)  specifies  that  the  field  amplitude  statistics  will  be  log-normal. 

Despite  the  success  of  the  Rytov  theory  and  related  mathematical  approximations,  there  are  many 
interesting  beam-propagation  scenarios  that  cannot  be  accurately  analyzed  in  this  manner.  The 
relationship  given  in  Eq.  (2-1 1)  is  immediately  suspect,  since  it  predicts  that  a  y  will  become  infinite  as 

L  — >  oo  .  This  result  is  inconsistent  with  simulations  of  large-diameter  transmitters,  which  show  that  cr y 
approaches  a  constant  as  the  propagation  range  becomes  large. 

Transmission  from  apertures  of  finite  dimension  is  known  as  beamwave  propagation,  and  for  many 
geometries  of  interest  the  associated  scintillation  properties  are  not  well  characterized  by  relationships 
derived  for  the  point-source  or  infinite  plane-wave  cases.  Furthermore,  numerical  simulations  and  field 
experiments  both  demonstrate  that  the  statistics  of  fading  due  to  propagation  through  strong  turbulence 
are  not  log-normal.  In  short,  many  of  the  assumptions  commonly  applied  to  the  design  of  systems  subject 
to  strong  channel  fluctuations  are  inappropriate  and  may  lead  to  serious  errors  in  performance  estimates. 
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3.  THE  THREE  STAGES  OF  OPTICAL  SCINTILLATION 


The  optical  field  relationships  given  in  Eqs.  (2-7)  and  (2-8)  suggest  a  multiplicative  process  that 
imposes  a  weak  distortion  onto  a  diffraction-limited  beam  profile.  Numerical  simulations  (which  are 
described  later  in  this  document)  verify  this  prediction  when  the  integral  of  the  Cn  profile  over  the 
propagation  path  is  small.  A  good  measure  of  the  absolute  turbulence  strength  is  the  coherence  diameter, 
/q,  defined  by  Fried  [9] 


/n  — 


0.423  kl 


\C2n{z)dz 


-1-3/5 


(3-1) 


It  can  be  shown  that  the  relative  impact  of  the  atmosphere  on  the  far-field  quality  of  a  beam  transmitted 
from  an  aperture  of  diameter  D  is  directly  related  to  the  D/r0  ratio.  This  quantity  plays  a  primary  role  in 
the  characterization  of  scintillation  behavior. 

Charnotskii  [10]  was  one  of  the  first  investigators  to  propose  a  segmented  classification  of 
scintillation  based  on  a  multidimensional  parameter  space,  and  his  article  suggests  that  as  many  as  1  I 
distinct  classes  are  needed  to  properly  describe  all  aspects  of  this  complex  phenomenon.  Charnotskii’s 
research  is  unique  in  that  it  concludes  that  a  thorough  analysis  of  scintillation  cannot  be  based  solely  on  a 
single  set  of  weak-turbulence  approximations.  One  would  hope,  however,  that  it  might  be  possible  to 
capture  the  essential  aspects  of  this  process  in  a  somewhat  simpler  set  of  mathematical  constructs. 

The  phenomenological  approach  proposed  in  this  document  segregates  scintillation  behavior  into 
three  classes,  and  associates  each  class  with  a  single  observable — the  D/r$  ratio.  This  methodology 

derives  from  a  detailed  analysis  of  vertical-path  and  horizontal-path  numerical  simulations  performed  at 
Lincoln  Laboratory.  On  a  qualitative  level  these  new  results  are  consistent  with  the  work  of  prior 
investigators,  and  the  predicted  numerical  values  for  very  weak  and  very  strong  turbulence  conditions 
were  correctly  reproduced.  (Contrary  to  the  prediction  that  obtains  from  Eq.  (2-11),  a %  achieves  a 

saturation  value  of  0.25  for  severe  channel  conditions.)  When  <j  x  is  plotted  as  a  function  of  Dfv§  ,  a  clear 

trend  is  revealed.  As  shown  in  Figure  3-1,  for  a  scenario  in  which  the  r0  value  is  fixed,  the  log-amplitude 
variance  is  constant  when  D/tq  «  1  and  approaches  the  saturation  limit  when  D/r#  » 1  .  The  transition 
region  occurs  when  Djr§  is  in  the  range  of  1  to  10.  This  behavior  has  been  observed  for  a  variety  of 
turbulence  profiles  and  for  both  vertical  and  horizontal-path  geometries. 
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Figure  3-1 .  Estimates  of  as  a  function  of  the  normalized  field  perturbation  strength  for  vertical  path 

propagation.  These  results  were  obtained  from  numerical  simulations  of  beam  wave  projection  through  a 
series  of  discrete  phase  screens  from  a  ground-level  transmitter  to  a  receiver  in  geostationary  orbit.  In 
the  small  D/Vq  limit  these  data  agree  with  the  Rytov  prediction  for  point-source  scintillation. 


While  plots  similar  to  that  shown  in  Figure  3-1  have  frequently  been  published,  little  attention  has 
been  paid  to  the  physical  processes  that  shape  the  scintillation  curve  beyond  the  region  where  the  Rytov 
approximation  is  applicable.  A  detailed  examination  of  far-field  beam  images  derived  from  numerical 
simulations  reveals  clear  distinctions  between  the  profiles  for  small,  intermediate,  and  large  values  of 
D/Vq  ;  these  differences  are  illustrated  by  the  examples  shown  in  Figure  3-2.  The  short-exposure  image 
shown  in  the  leftmost  picture  indicates  the  behavior  when  D/r$  is  small.  In  this  case  the  long-exposure 
beam  width  is  approximately  diffraction-limited,  but  individual  realizations  display  a  relatively  weak 
multiplicative  noise  component.  This  is  the  only  region  in  which  the  Rytov  approximation  is  strictly  valid 
and  for  which  the  distribution  of  intensity  fluctuations  is  log-normal. 


% 
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D/r0  «  1  Z)/r0  «  1  £>/r0  »  1 


Figure  3-2 .  Beam  profile  evolution  as  a  function  of  the  normalized  field  perturbation  strength.  These 
images  indicate  three  distinct  physical  mechanisms  for  irradiance  fluctuations  in  the  plane  of  the 
receiver.  The  tilt-jitter  effects  arid  beam  breakup  shown  in  the  second  and  third  images  are  not  accounted 
for  by  the  Rytov  theory . 


When  D/r0  is  of  order  unity,  the  strength  of  the  log-normal  fluctuations  begins  to  diminish,  and  the 
scale  size  of  the  noise  peaks  becomes  comparable  to  the  diffraction-limited  beam  diameter.  However,  tilt 
jitter  becomes  more  pronounced,  and  low  spatial-frequency  aberrations  distort  the  short-exposure  beam 
shape.  Both  of  these  effects  result  in  a  dramatic  increase  in  far-field  irradiance  fluctuations.  Although 
gross  beam  motion  is  not  usually  thought  of  as  scintillation,  the  distinction  is  artificial  since  both 
introduce  a  corruption  of  the  on-target  irradiance.  The  statistical  properties  for  this  case  are  not  log¬ 
normal  and  are  strongly  dependent  on  the  quality  of  the  track  correction  used  in  the  transmitter. 

Tilt  fluctuations  are  associated  with  the  appearance  of  the  lowest  spatial-frequency  phase-distortion 
components,  which  are  the  dominant  distortion  factor  [11].  As  the  strength  of  the  higher  spatial-frequency 
constituents  becomes  comparable  to  2 k  ,  the  beam  sidelobes  increase  in  intensity  until  the  mainlobe  and 
the  sidelobes  become  indistinguishable.  Full  saturation  occurs  when  the  phase  errors  across  the  aperture 
are  uniformly  distributed  between  -k  and  n  .  Negative-exponential  statistics  are  observed  for  these 
conditions. 

The  statistical  descriptions  developed  in  the  following  sections  are  a  direct  product  of  observations 
drawn  from  Figures  3-1  and  3-2.  The  intent  of  this  effort  is  to  develop  a  better  understanding  of  the  causes 
and  the  effects  of  atmospheric  turbulence  beyond  the  range  in  which  the  Rytov  approximation  is  valid. 
The  goal  is  to  provide  a  set  of  mathematical  constructs  that  can  be  used  to  predict  system  performance  for 
any  propagation  geometry  and  turbulence  environment. 
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4.  A  NEW  PERSPECTIVE  ON  THE  BEAMWAVE  SCINTILLATION  PROBLEM 


The  most  important  observation  made  in  the  previous  section  is  that  the  Rytov  approximation  fails 
for  the  scenarios  of  greatest  concern  to  the  designers  of  free-space  optical  systems  that  are  intended  to 
function  in  a  stressing  turbulence  environment.  Even  a  cursory  analysis  of  numerical  results  for  moderate 
to  strong  perturbations  reveals  serious  discrepancies  between  predictions  and  measurements.  To  develop 
improved  performance  models,  the  effects  of  tilt  jitter  and  beam  breakup  must  be  explicitly  incorporated. 

In  the  investigation  of  this  problem,  it  was  soon  realized  that  diffraction  plays  a  key  role  in 
beamwave  scintillation.  The  location  of  the  perturbation  field  relative  to  the  Rayleigh  range,  where 
diffractive  beam  spreading  becomes  important,  strongly  impacts  the  strength  and  statistical  behavior  of 
signal  fluctuations  measured  at  the  receiver  plane.  It  is  common  to  refer  to  field  points  that  are  within  the 
Rayleigh  range,  =  D2// 1 ,  as  being  in  the  near  field  of  the  transmitter,  and  describe  more  distant  points 
as  lying  in  the  far  field.  For  the  case  of  beam  propagation  from  a  large  transmitter  from  ground  to  space, 
turbulence-induced  distortions  will  be  solely  confined  to  the  near  field  of  the  source.  As  the  transmitter 
diameter  shrinks,  the  Rayleigh  range  moves  closer  to  the  source  and  the  accumulated  perturbations 
migrate  to  the  far  field.  In  the  first  example,  the  fields  in  the  receiver  and  transmitter  planes  form  a  simple 
Fourier-conjugate  pair.  A  quantitative  description  of  the  latter  case  is  significantly  more  complex,  and  this 
difference  is  central  to  the  discussions  that  follow. 

Although  this  analysis  is  intended  to  be  as  general  as  possible,  three  simplifying  constraints  will  be 
applied  that  are  applicable  to  most  system  constructs  of  practical  utility.  The  first  is  that  the  receiver 
aperture  is  small  compared  to  the  size  of  the  speckle  features  in  the  far-field  beam  and  that  it  lies  along 
the  nominal  propagation  path.  This  assumption  enables  scintillation  analyses  to  be  performed  on  a  single 
point  within  the  plane  of  the  receiver.  The  second  constraint  is  that  the  receiver  lies  at  the  focus  of  the 
transmitted  beam.  For  most  systems  of  practical  interest,  the  average  received  irradiance  is  maximized 
when  the  focus  applied  to  the  transmitted  beam  is  equal  to  the  receiver  range;  in  addition,  this  geometry 
also  establishes  an  exact  Fourier-transform  relationship  between  the  two  optical  planes.  The  third 
assumption  is  that  the  transmitted  beam  is  an  infinite  Gaussian  characterized  by  the  field  radius,  wp. 
Although  much  of  the  turbulence  literature  assumes  transmission  from  a  uniformly-illuminated  aperture 
of  diameter  D,  the  beams  from  most  transceivers  are  best  described  by  a  truncated-Gaussian  shape.  It  will 
be  shown  that  relationships  developed  for  uniformly-illuminated  apertures  can  be  easily  converted  to 
Gaussian-beam  expressions  through  the  application  of  the  scaling  transformation,  D  -  2 wq  . 

4.1  BEAMWAVE  PROPAGATION  RELATIONSHIPS  FOR  UNPERTURBED  GAUSSIAN 
BEAMS 

The  evolution  of  the  beamwave  profile  with  propagation  range  is  shown  graphically  in  Figure  4-1. 
The  spatial  and  temporal  characteristics  of  irradiance  fluctuations  resulting  from  perturbations  introduced 
in  the  near-  and  far-field  regions  of  the  beam  differ  significantly,  and  it  is  necessary  to  develop  a  separate 
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set  of  analytical  tools  for  each  region.  The  basic  mathematical  relationships  that  are  later  used  to  construct 
statistical  irradiance  models  are  developed  in  this  section. 


Beamwave 


< - 

Transmitter 

Plane 


L 


Receiver 

Plane 


Figure  4-1.  Evolution  of  a  focused-beam  profile  from  transmitter  to  receiver  planes .  The  shape  of  the 
beam  is  governed  by  geometrical  propagation  in  the  near-field  of  the  transmitter  and  by  diffractive 
spreading  in  the  far  field .  The  transition  between  these  two  regions  occurs  at  the  Rayleigh  range. 


The  composite  beamwave  profile  shown  in  Figure  4-1  shows  that  the  shape  of  the  focused  beam  is 
established  by  the  geometrical  propagation  pattern  in  the  transmitter’s  near-field  region,  but  becomes 
diffraction-limited  in  the  far-field.  A  Gaussian  beam  of  infinite  extent  remains  Gaussian  over  the  entire 
range  and  is  completely  characterized  by  the  range  dependent  radius,  vv(r) ,  and  wavefront  focal  radius, 
f(r).  The  following  derivation  borrows  heavily  from  Goodman’s  discussion  of  Fresnel  propagation 
[12]. 


The  expression  for  the  optical  field  amplitude  in  the  plane  perpendicular  to  the  propagation 
direction  derives  from  the  Huygens-Fresnel  principle 


OO  00 

u\(x\>y\)=  J  |  h(x0,y0-,x],y\)U0(xQ,y0)dxQdy0  ,  (4-1) 

— oo  —oo 
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where  U0  is  the  complex  field  amplitude  within  the  x,y  plane  at  range  z(h  and  Ui  is  the  amplitude  at  z/.  The 
propagation  transfer  function,  h ,  is  given  by 


ik 


2/r  (z\  -Zq) 


exp  ji  k  yj(x\  -  x0)2  +  (.v,  -  vq)2  +  (z,  -  z0 )2  j  , 


(4-2) 


3  7T  2  2 

where  k-ln/X.  For  most  problems  of  interest  (z j  -zq)'  » —  (x(  -xq)  +  (^|  -j'q) 

4/1  L 


,  in  which 


case  the  Fresnel  approximation  is  applicable 

i  k  exp|/ k  (z|  -  z0) 


M-WoiWi)38-' 


exp  <:  i  - 


2^(z,-z0)  [  2(z, -z0)L 


(x,  -x0)2  +  (y,  -y0f 


(4-3) 


An  important  property  of  the  Fresnel  approximation  is  that  in  the  absence  of  field  perturbations  a 
single  propagation  step  from  zo  to  z 2  yields  the  same  result  as  multiple  steps  within  the  same  interval. 
Consider  the  two-step  process  zq  =>  z j  =>  zj  described  by  the  equation 
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Observe  that 
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where  a  =— - -  and  b  =  — - and  note  that  the  integral  over  x}  can  be  performed  using  the 

2{z\  ~z0)  2{z2-z\) 

identity 
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This  yields  the  result 
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With  the  exception  of  the  leading  phase  constant,  this  is  equivalent  to  the  propagation  expression  for  the 
single-step  propagation  zq  — >  z2  . 

jt  i  2  2  \ 

When  the  propagation  range  is  extended  further  such  that  Z\  -z0  »  —  Lvq  +  v0  ,  the 

X  '  /max 

Fraunhofer  approximation  provides  an  accurate  estimate  of  the  transformed  field 
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(4-4) 


Note  that  aside  from  phase  terms  that  are  independent  of  x{)  and  y(>,  this  expression  is  equivalent  to  a 
Fourier  transform.  Goodman  shows  that  this  approximation  is  also  valid  for  short  propagation  distances  if 
the  transmitted  beam  is  focused  at  z{.  It  is  also  worth  noting  that  along  the  boresight  path  the  field 
expression  for  Fraunhofer  propagation  is  proportional  to  the  two-dimensional  integral  of  U(l 


U\  (0,0) = 


2^-(z,-z0)  J 


J  J  uo{xo>yo)dxodyo 


(4-5) 


which  is  a  useful  relationship  for  numerical  calculations. 

The  field  transformation  concepts  just  presented  provide  a  basis  for  a  discussion  of  the  effects  of 
turbulence-induced  phase  perturbations.  An  important  element  of  the  following  analysis  is  the 
demarcation  between  the  near  field  and  far  field  of  the  transmitted  beam  and  the  placement  of  the 
perturbation  with  respect  to  that  transition  point.  The  beamwave  problem  can  be  significantly  simplified 
by  restricting  the  discussion  to  focused  Gaussian  beams.  Investigations  by  Andrews  and  Phillips  have  also 
concentrated  on  Gaussian-beam  behavior  [13]. 

For  virtually  all  cases  of  interest,  optimal  performance  is  achieved  when  the  geometrical  focal 
distance  and  the  range  to  the  receiver  are  identical.  The  field  at  the  transmitter  will  be  defined  as  follows 

(ro )  -  4)  exp {-*'  k  fQ  ftfo }  exp {-r02 / wq  }  ,  (4-6) 

where  wn  is  the  field  radius  and  f{}  is  the  focal  distance.  Assuming  Fresnel  propagation  from  z^O  to  z/,  the 
field  expression  at  range  z/  is 
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The  integral  relationship  exp  j-«  x2 / 2^  Jq(x y)  dx  =  —  exp  j -y2 /la 


can  be  applied  to  Eq.  (4-7),  where  x  =  r$ ,  y  =  k  >\  /z, ,  and  #= — — -/A 
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This  field  description  has  the  same  structure  as  the  original  beam,  in  that  both  the  amplitude  and  phase 
profiles  are  Gaussian.  The  beam  radius  at  range  z}  is 


4zi2  +  k2 
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(4-9) 


where  the  first  term  is  the  diffractive  beam-spreading  factor  and  the  second  describes  the  beam 
dimensions  for  geometrical  propagation.  The  wavefront  is  defined  by  the  focal  distance 


f\  =  z! 
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(4-10) 


Note  that  the  beam  radius  at  range  z/  is  minimized  when  /0  =  z\ . 
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It  will  be  useful  to  define  a  range  parameter  that  specifies  the  onset  of  beam  diffraction.  According 
to  Eq.  (4-9),  the  beam  radius  will  increase  by  a  factor  of  \fl  when  the  propagation  range  is 


ZR  “ 
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2 


,  k  wq 
2/o 


(4-11) 


To  compare  this  result  with  derivations  developed  for  uniformly-illuminated  apertures  of  diameter,  D,  it  is 
useful  to  establish  a  relationship  between  D  and  w(}.  When  the  Airy  irradiance  pattern  derived  from  the 
Fourier  transformation  of  a  circular  aperture  is  compared  with  the  far-fleld  profile  of  a  diffraction-limited 
Gaussian  beam,  the  main-lobe  structures  are  found  to  be  closely  matched  when 

D  =  2i/I  w0  .  (4-12) 


Therefore,  fora  uniformly  illuminated  circular  aperture 
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For  /q  =  oo  this  result  closely  approximates  the  standard  definition,  zR  =  D" / A  . 
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Well  beyond  the  Rayleigh  range  the  beam  assumes  the  characteristics  of  a  spherical  wave 
originating  from  a  point-source  transmitter,  as  illustrated  earlier  in  Figure  4-1.  In  the  near  field,  where 
Z|  «  zr  ,  the  approximate  expressions  for  w/,//,  and  the  multiplicative  constant,/!/,  are 
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whereas  in  the  far  field  of  the  transmitter 
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Once  the  far-field  range  limit  has  been  reached,  the  beam  radius  and  wavefront  focus  both  scale  linearly 
with  propagation  distance. 
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4.2  RECEIVER-PLANE  FIELD  EFFECTS  DUE  TO  NEAR-FIELD  PHASE  PERTURBATIONS 


For  a  few  propagation  geometries,  such  as  vertical  link  between  a  large  ground-based  telescope  and 
a  high-altitude  satellite,  the  accumulated  channel  distortions  all  effectively  lie  in  the  plane  of  the 
transmitter.  In  this  case,  the  phase  errors  for  a  single  realization  can  be  added  linearly  and  treated  as  a 
single  distortion  screen  defined  by  the  function,  ^(x0,yo)-  The  transmitter-plane  representation  of  a 

perturbed  Gaussian  beam  focused  on  a  receiver  at  range  z-L  is 

O0(x0,  To )  =  A)  exP {*  4> (*0  ’  To )}  exP \-ik  (xo  +  To  )I2L)  exP  f -  (-*0  +  To  )/**£ }  •  (4-1 6) 


The  inclusion  of  the  focus  term  places  the  receiver  in  the  Fourier  transform  plane  of  the  transmitter,  in 
which  case  the  Fraunhofer  approximation  can  be  applied  to  obtain  the  field  at  range  L 
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(Complex  constants  having  unity  amplitude  have  been  eliminated  from  this  expression  to  simplify  the 
notation.)  With  the  variable  substitutions,  a  =  kx^/L  and  -k  yo/L  ,  the  previous  equation  can  be 

written  in  the  form  of  a  Fourier  transform 


0L  {xL,yL )  =  expj/J:  (xj  +  y\ )/2uJ 
x  y  lexpj i</>(La/k,L/3/k)j  exp 


L2a2+L2y)/k2^  j}. 


(4-18) 


which  explicitly  shows  the  range  and  wavelength  dependence  of  the  radial  scaling.  Note  that  the  leading 
phase  term  is  equivalent  to  that  for  a  point  source  located  in  the  plane  of  the  transmitter,  and  in  the 
absence  of  phase  errors  the  second  term  is  a  simple  Gaussian  amplitude  profile  as  described  previously  in 
Eq.  (4-8).  It  should  also  be  emphasized  that  since  the  outgoing  wavefront  has  been  tailored  to  optimize 
received  power,  all  phase  distortions  inserted  in  the  near-field  result  in  a  degradation  of  the  far-field 
signal. 

Recall  that  this  derivation  is  based  on  the  assumption  that  all  channel  effects  are  confined  to  the 
plane  of  the  transmitter.  This  constraint  can  be  relaxed  somewhat  by  allowing  the  field  perturbations  to  be 
placed  further  downrange,  but  within  the  Rayleigh-range  limit.  If  the  screen  is  located  at  range,  z,  with 
respect  to  the  transmission  aperture,  then 


17 


(4-19) 
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where  y  represents  the  radial  scaling  relationship  between  the  unperturbed  beam  at  the  transmitter  and  the 
unperturbed  beam  at  range  z.  Thus,  from  Eq.  (4-9) 
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Since  the  beam  is  focused  on  the  receiver,  it  initially  shrinks  in  dimension  but  increases  in  size  as  the 
Rayleigh  range  is  approached.  The  assumptions  of  the  near-field  approximation  are  clearly  invalid  when 
diffraction  effects  start  to  dominate  and  ^becomes  significantly  larger  than  1.  On  the  other  hand,  one 
could  treat  the  previous  expression  as  the  sum  of  a  geometrical-propagation  component,  which  describes 
a  beam  that  is  focused  to  a  point  at  the  receiver  and  a  spherical-wave  component  that  emanates  from  the 
transmitter.  The  first  of  these  two  terms  yields 


r * 


(4-2  l ) 


which  has  appeared  in  a  number  of  prior  publications  in  which  scintillation  effects  for  focused-beam 
geometries  have  been  investigated  [8].  An  accurate  treatment  of  the  spherical-wave  component  appears  to 
require  a  very  different  set  of  assumptions  and  approximations,  as  described  in  the  following  section. 

4.3  RECEIVER-PLANE  FIELD  EFFECTS  DUE  TO  FAR-FIELD  PHASE  PERTURBATIONS 


When  the  diameter  of  the  transmitter  becomes  small  enough  to  push  most  of  the  channel  effects 
beyond  the  Rayleigh  range,  the  propagation  geometry  begins  to  approach  the  point-source  limit.  In  this 
case  the  perturbation  and  receiver  planes  cannot  be  treated  as  a  Fourier  conjugate  pair,  and  their 
relationship  must  instead  be  defined  through  the  Fresnel  transformation.  In  general,  signal  fluctuations 
due  to  far-field  perturbations  are  less  severe  since  the  relative  separation  between  the  two  planes  is 
smaller.  A  second  important  difference  between  the  near-field  and  far-field  geometries  is  that  in  the  latter 
case  the  insertion  of  a  phase  screen  will  result  in  localized  regions  of  signal  enhancement  due  to  focus 
contributions  that  counteract  diffractive  wavefront  expansion  (see  Figure  4-2);  when  the  channel  effects 
are  sufficiently  weak,  the  mean  irradiance  will  be  identical  to  the  diffraction-limited  value.  Both  of  these 
factors  make  it  easier  to  apply  Rytov  analysis  to  far-field  distortions. 
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Plane 

Figure  4-2 .  Illustration  of  the  effect  of  placing  a  focusing  element  in  the  fa r-field  of  a  point  source.  A 
point  along  the  propagation  axis  will  experience  an  increase  or  decrease  in  the  optical  irradiance, 
depending  on  the  sign  of  the  focusing  element.  Note  that  phase  distortions  that  result  in  an  on-axis 
irradiance  enhancement  only  occur  in  the  far  field  of  the  transmitter. 


A  mathematical  description  of  Fraunhofer  propagation  from  a  source  at  z-  0,  through  an 
intermediate  plane  at  Zj,  and  then  to  a  second  optical  plane  at  z2  vvas  provided  in  Section  4.1.  To 
simplify  the  notation  and  maintain  consistency  with  the  previous  discussion,  in  this  section  the  receiver’s 
range  will  be  represented  as  L,  and  the  range  to  the  intermediate  perturbation  plane  will  be  given  as  z. 

If  the  assumption  is  made  that  z  »  ,  a  Gaussian  beam  that  encounters  the  perturbation  field,  fa 

can  be  described  as  follows 


0  (x,  y)  =  A  exp|/^(x,y)j  exp  j//:|x2  +  v2  j^2zj  exp  j-|x2  +  y2 jy/w2  |  , 


2  z  k  n 

and  as  specified  in  Eq.  (4-15)  ww -  ,  and  A  ~ 


k  w, 


'0 


2  z 


(4-22) 

(4-23) 
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Fresnel  propagation  to  L  yields 


^l{xl^l)  = 


2n  (L  -  z) 


QO  oc 

r  r 

Lrb-^)h(frv)2 

J  J  exP< 

-00  — oc- 

[  ‘  2(i-0  j 

U(x,y)dxdy  .  (4-24) 


This  expression  can  be  rearranged  to  incorporate  the  focus  term  in  the  Fresnel  kernel,  so  that 

k 


Ul{xL^l)  = 


2tt  (L -z) 


exp  \i  k  \  x\  +  y2L  1/2  L 


OC  00 

V'  1  p  v  r"\  J 

\  L  (zxjL-xf  +(zyjL-yf 

x  exp  > 

—00  —00 

[  z  2  {L-z) 

(4-25) 


£/'(jc,  y)  dx  dy 


where 


Uf(x,y )  =  zlexpj/ <^(x,  exp  j-|x2  +  y 


2  2 
W 


(4-26) 


is  the  product  of  the  Gaussian  amplitude  and  the  perturbation  field  atz.  The  field  transfer  function  can  be 
expressed  as 


h'(x,y\xL,yL) 


2  x(L~; 


-exp- 


(  2  ,  2)1 

r 

[XL  +);L  j 

>exp-j 

2  L 

J 

L  {zxJ_L-xf  +{zyj_  ,'L- _vf 
z  2  (L-z) 


(4-27) 

In  the  limit  that  the  propagation  distance  L-z  is  small  compared  to  k  w  ,  the  third  term  of  the  previous 
expression  can  be  replaced  by  a  pair  of  delta  functions  by  applying  the  identity 

£(*)slim  exp \-ax2\  Jrtfa  , 

a— >0  v  n 

thus  yielding  the  result 


h'(x,  y\xL,yL)  w^-exp 


i  k  [x\ \  +  y\ ) 


2  L 


>  S\  — x,  -  X 
L  L 


J 


■yL 


(4-28) 


To  first  order,  the  field  at  L  is  simply  a  replica  of  (/,  which  has  been  scaled  along  the  radial  axes  by  the 
factor  Ljz 
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(4-29) 


0L(xL,yL)*j-<*?\ik[xl+yl)/2L\W 


This  solution  represents  the  trivial  case,  where  the  beam  propagation  distance  beyond  the  phase  screen  is 
too  short  to  substantially  modify  the  receiver-plane  irradiance  profile. 

As  the  separation  distance  between  z  and  L  increases,  the  delta  functions  that  appear  in  the  previous 
expression  must  be  replaced  with  a  more  accurate  approximation.  As  illustrated  in  Figure  4-3,  the 

integrals  of  both  sin^aarj  and  cos|ajc2j  display  an  abrupt,  and  approximately  linear,  transition 

centered  at  x  =  0  .  The  derivatives  of  both  functions  can  therefore  be  replaced  by  the  rectangle  function, 
Yliyjlaj K  Jcj ,  where  the  normalization  is  given  by  the  integral 


OO  CO 

exp j/a  (T  +  y2  jj  dbcJy  =  /  —  . 

—00  -OO  ^ 

The  rectangle  function  establishes  the  linear  dimensions  of  a  characteristic  integration  region 


l  =  . 


2 n  z{L-: 

k  1 


(4-30) 


This  results  in  the  replacement  of  the  delta  functions  with  a  two-dimensional  average  of  U'(x,y) 
centered  at  the  point  [z  x^/ L,zy^/ L\ ,  i.e., 


0L  (xl  >  -vl  )  *  A  j  exP  {ik(xl+  yl)/2L] 

1  i/2  1/2  \  ( z  z  X 

J  |  exp\i<f>\-xL-x,-yL-y 


i 


L 


L  '  ) 


-exp 


w 


L 


xL-x\  + 


-jyL-y 

VT  y 


-t/2-t/2 

At  the  nominal  boresight  location,  the  expression  is  further  simplified 

111  1/  2 

0L( 0,0)»-^-  J  J  exp|/^(x:,_p)j  expj-^x2  +  t2)//vv2|  dxdy  , 
Z  L  _e/2-£/2 

and  if  it  is  assumed  that  the  physical  dimensions  of  the  beam  at  z  are  much  larger  than  £ 


(4-31) 


>dxdv . 


(4-32) 
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(4-33) 


t/2  £/2 

j  {  exp {i((>{x,y)}dxdy 
^  L  -t/2-  t/2 

For  the  statistical  analysis  developed  in  the  next  section,  it  will  be  convenient  to  replace  £  with  the 
circular  diameter 


or  equivalently  the  Gaussian  radius 


we  -  2  ^ 2  dl  - 


(4-34) 


The  physical  interpretation  of  this  result  is  that  the  scintillation  experienced  by  a  point  in  the 
receiver  plane  is  only  sensitive  to  far-field  channel  distortions  that  lie  along  the  propagation  axis  and  are 
within  a  radius  defined  by  Wg  .  However,  as  the  range  to  the  perturbation  plane  approaches  the  near-field 
region,  this  simple  analytical  description  must  be  modified  to  account  for  plane-wave  beam  propagation 
within  the  transmitter  radius,  Wq.  As  discussed  below,  these  factors  significantly  reduce  the  impact  of 
near-field  phase  fluctuations  on  the  far-field  contribution  to  the  net  scintillation. 


Figure  4-3.  Comparison  of  the  integrals  of  cos  (xf  and  sinfx")  with  the  integral  of  the  rectangle  function. 
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4.4  NUMERICAL  SIMULATION  TOOLS 


The  discussion  of  scintillation  statistics  provided  in  the  next  section  derives  from  the  analytical 
analysis  presented  earlier  in  this  section,  but  is  also  heavily  dependent  on  the  results  of  numerical 
simulations  of  beamwave  propagation  performed  at  Lincoln  Laboratory.  The  current  version  of  the  code, 
known  as  POPS,  is  a  descendent  of  the  three-dimensional  time-dependent  program,  MOLLY,  which  was 
written  at  Lincoln  Laboratory  in  the  late  1980’s  to  study  turbulence  and  thermal  blooming  effects.  Both 
programs  incorporate  subroutines  to  generate  random  Kolmogorov  phase  screens  [14]  and  a  means  to 
propagate  from  one  screen  to  the  next  using  a  two-dimensional  Fresnel  transformation  [12].  As  the  beam 
propagates  from  screen  to  screen,  the  Talanov  scale  transformation  [15]  is  imposed  to  maintain  a  constant 
ratio  between  the  beam  size  and  the  dimensions  of  the  field  discretization  grid.  An  arbitrary  number  of 
phase  screens  can  be  created  to  simulate  the  desired  Cn  profile,  and  each  screen  can  be  moved 

perpendicular  to  the  beam  in  accordance  with  the  wind  velocity  profile  for  time-dependent  simulations. 
Statistical  behavior  is  typically  characterized  by  ainning  100,000  independent  realizations  of  the 
propagation  channel.  POPS  is  a  state-of-the-art  C-based  code  with  a  web-browser  interface,  which 
currently  functions  on  a  16-node  Xeon  processor  cluster. 

For  statistical  investigations,  the  output  of  the  POPS  code  is  generally  configured  to  record  the 
irradiance  at  a  small  number  of  grid  points  within  the  receiver  plane.  However,  to  obtain  a  better 
understanding  of  the  physical  processes  that  drive  beam  scintillation,  it  is  often  useful  to  study  selected 
images  of  the  entire  profile.  The  evolution  of  the  dominant  phenomenological  mechanisms  that  drive  far- 
field  scintillation  are  clearly  illustrated  by  the  beam  images  shown  in  the  next  section.  In  support  of  this 
investigation,  an  extensive  set  of  simulations  have  been  performed  for  the  vertical  propagation  geometry 
to  a  receiver  at  40  Mm  and  horizontal  propagation  to  a  receiver  at  10  km.  The  tilt-removed  results  that  are 
described  in  the  next  section  were  obtained  by  measuring  the  irradiance  for  each  realization  at  the  grid 
point  located  at  the  center  of  mass  of  the  far- field  beam  profile. 
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5.  IRRADIANCE  STATISTICS  AT  THE  RECEIVER-PLANE  BORESIGHT 


In  this  section  separate  scintillation  index  expressions  will  be  developed  for  each  of  the  three  beam- 
propagation  regimes  defined  previously;  specifically,  w0/r0  «1,  w0/r0  *] ,  and  wq/Vq  »  1  .  It  is  only 

in  the  first  of  these  cases,  where  the  geometry  approaches  the  point-source  limit,  that  far-field  channel 
effects  tend  to  dominate  scintillation  behavior. 

The  parameter  terminology  introduced  in  the  previous  section  will  be  continued  in  this  discussion, 
including  the  use  of  the  Gaussian  beam  radius  Wq  to  characterize  the  dimensions  of  the  transmitted  beam. 
This  analysis  begins  with  the  case  that  is  best  characterized  by  traditional  statistical  descriptions  based  on 
the  Rytov  approximation. 

5.1  SCINTILLATION  INDEX  FOR  w0/r0  «  1 :  WEAK  FAR-FIELD  TURBULENCE 

The  propagation  relationships  developed  in  Section  4.3  demonstrate  that  when  the  transmitted  beam 
diameter  is  very  small,  the  phase  perturbation  and  the  receiver  planes  both  lie  in  the  far  field  of  the 
transmitter.  In  this  limit,  the  detected  irradiance  fluctuations  due  to  a  perturbation  field  at  range  2  are  only 
influenced  by  an  illuminated  region  of  radius 


Wo 


R~z) 

kL 


(5-1) 


which  is  generally  much  smaller  than  the  physical  dimensions  of  the  diffraction-limited  beam.  To  first 
order,  the  net  impact  of  phase  distortions  within  this  region  can  be  modeled  as  a  thin  lens  that  either 
enhances  or  counteracts  diffractive  beam  spreading.  At  the  receiver  plane,  this  effect  produces  a  weak 
multiplicative  modulation  of  the  diffraction-limited  beam  profile,  as  quantified  earlier  in  Eq.  (4-33). 

The  beam  image  and  profile  examples  given  in  Figure  5-1  show  little  evidence  of  either  pointing 
jitter  or  time-averaged  shape  distortions.  The  primary  feature  illustrated  in  these  pictures  is  a  high  spatial- 
frequency  speckle  pattern  that  imposes  relatively  shallow  fades.  In  this  environment,  the  assumptions 
associated  with  the  Rytov  approximation  are  valid,  and  numerical  simulations  show  that  the  on-axis 
irradiance  distribution  is  approximately  log-normal  with  a  mean  Strehl  that  is  close  to  unity. 
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RADIAL  DISPLACEMENT  (km) 


Figure  5-1.  Far-field  beam  image  (left  figure)  and  horizontal-axis  profile  (right  figure)  for  Wq/zq  =  0.01 

The  Hufnagel-  Valley  turbulence  profile  was  used  in  this  simulation  of  vertical  propagation,  which 
resulted  in  an  r()  value  of  19  cm  at  1 .55  mm.  The  propagation  distance  is  40  Mm,  and  the  linear  dimension 
of  the  image  fields  is  50  km.  The  long-exposure  beam  profile  is  approximately  diffraction-limited. 


Estimates  of  the  scintillation  index  for  the  point-source  limit  follow  directly  from  the  work  of  a 
number  of  investigators  in  the  field  of  adaptive  optics  [16].  Following  Sasiela’s  derivation  [8],  the 
focused-beam  scintillation  index  for  Kolmogorov  turbulence  would  be  expressed  as  follows 


aj  =0.829  k2 


k-"»  sin2 


k2  z{L-z) 
2k  L 


cIk  dz 


(5-2) 


where  k  is  the  spatial  frequency  vector.  The  sin 2  [•]  term  that  appears  in  this  expression  represents  the 
filter  function  for  a  spherical  aperture,  which  captures  the  effect  of  the  far-field  integration  region,  w{  , 
defined  in  Eq.  (5-1).  For  a  point  source,  the  evaluation  of  the  integral  over  the  zc  vector  gives 

aj  =  2.25  k7/6  jcj  (z)  z5/6  f— -1  dz  .  (5-3) 

o  \  L  J 


Several  authors  define  the  propagation  parameter,  y  =  zj L  ,  and  express  this  equation  as  follows 
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(5-4) 


a]  =22ik1'i)c1a(z)(L-zf6  /'<•  dz 

0 


The  previous  result  can  be  extended  to  include  transmitters  of  finite  dimension  by  appropriately 
modifying  the  argument  of  the  integral  in  Eq.  (5-2).  The  approach  employed  here  is  similar  to  that 
developed  by  Sasiela  to  evaluate  the  magnitude  of  scintillation  for  extended  sources  and  receivers.  The 
effects  of  an  extended  source  of  radius,  w’5 ,  situated  at  range  z  are  described  by  the  filter  function 


F(K)  = 


(5-5) 


so  that  Eq.  (5-2)  can  be  rewritten  as 

o]  =  0.829  k2  Jc^(z)  | at-1  ,//3  sin2 
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(Ik  dz  .  (5-6) 


When  the  ratio,  k  wf / z  ,  is  large,  the  integration  over  k  can  be  approximated  by 
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erf  *1.53 


K2  (9 


f  L  ^ 
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W, 


?/3  dz  . 


(5-7) 


The  transition  between  Eqs.  (5-3)  and  (5-7)  can  be  accomplished  by  postulating  a  range-dependent  source 
radius  of  the  form 


W's,  =  Wq 


1 1  +  0.719 


_L _ z_ 

{L-z)  k  Wq 


(5-8) 


which  is  proportional  to  Vz  beyond  the  Rayleigh  range.  The  final  result  can  be  expressed  in  the  form  of 
Eq.  (5-4),  but  with  the  following  modification  of  the  propagation  parameter 


r  = 


1  +  1.39 


(L-A 
v  L  ) 


(5-9) 


Note  that  with  this  version  of  the  y  parameter,  the  scintillation  index  approaches  zero  monotonically  as 
wq  becomes  large.  It  must  be  emphasized,  however,  that  this  prediction  assumes  that  all  perturbation 
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effects  are  confined  to  the  far-field  region  of  the  transmitter.  This  relationship  is  similar  to  an  expression 
derived  by  Sasiela  for  beamwave  propagation  [8].  An  alternative  derivation  of  the  scintillation  index  for 
this  case  is  given  in  the  Appendix. 


For  future  studies  of  time-dependent  behavior  it  is  perhaps  of  interest  to  briefly  address  the 
characteristic  feature  dimensions  of  the  speckle  pattern  illustrated  in  Figure  5-1.  An  examination  of 
sample  profiles  for  a  range  of  Wq  values  reveals  that  the  speckle  dimensions  are  relatively  insensitive  to 

the  size  of  the  transmitted  beam,  so  that  the  number  of  speckles  within  the  diffraction-limited  far- fie  Id 
profile  is  inversely  proportional  to  wq  .  This  result  is  supported  by  a  calculation  of  the  coherence  diameter 
measured  in  the  plane  of  the  receiver  for  small-beam  illumination.  Using  the  definition  for  /given  above, 
the  characteristic  angular  dimension  of  the  speckle  features  is  predicted  to  be 


0.423  k2  |  C2 (z)  z5/3 
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1  +  1.39 
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Wn 


V-7/3 


(5-10) 


For  small  values  of  wq  this  result  is  independent  of  the  transmitted  beam  diameter,  and  the  estimates 
derived  from  this  expression  are  found  to  be  in  good  agreement  with  profile  realizations  obtained  from 
numerical  simulations. 


Figure  5-2  shows  a  plot  of  scintillation  index  values  derived  from  the  POPS  code  for  an  uplink 
geometry  through  a  Flufnagel-Valley  turbulence  profile  and  horizontal-path  propagation  over  a  10  km 
path  for  C2  =  1 .39x1 0~ 1 6  m-2  3  .  In  both  cases,  the  r0  value  was  computed  for  focused-beam  propagation 
from  0  to  L,  as  described  in  Eq.  (3-1) 


rC\  = 


0.423  A:2  fc2(z) 


L-z 


x5/3 


1  -V5 


dz 


\  z  ) 


The  dashed  lines  included  in  this  figure  indicate  the  scintillation  index  predicted  for  the  small  wq/tq  limit 

using  the  propagation  parameter  expression  given  in  Eq.  (5-9).  The  model  accurately  describes  the  data 
until  wq/Vq  approaches  unity. 

All  of  the  results  presented  in  this  section  are  consistent  with  traditional  analyses  based  on  the  first- 
order  Rytov  approximation.  In  the  next  two  sections,  it  will  be  shown  that  alternative  methodologies  are 
required  when  wq/tq  becomes  large. 
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Figure  5-2.  Comparison  of  simulation  results  (red  and  blue  data  points)  with  the  model  for  turbulence 
concentrated  in  the  far-jield  of  the  transmitter  (blue  dashed  line).  The  upper  figure  represents  an  uplink 
geometry  and  the  lower  is  horizontal  propagation  to  10  km. 
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5.2  SCINTILLATION  INDEX  FOR  w0/r0  »  1  :  MODERATE  NEAR-FIELD  TLRBULENCE 


In  many  respects,  the  moderate  turbulence  case  is  the  most  interesting  and  least  understood  of  the 
three  scenarios  treated  in  this  report.  Numerical  simulations  performed  at  Lincoln  Laboratory  show 
significant  discrepancies  between  measured  statistical  behavior  and  the  predictions  based  on  log-normal 
models.  A  distinguishing  characteristic  of  the  irradiance  distributions  derived  from  large  ensembles  of 
independent  channel  realizations  is  an  extended  power-law  tail  in  the  low-signal  region  of  the  probability 
density  function.  The  impact  of  this  discrepancy  is  that  the  log-normal  models  grossly  underestimate  the 
likelihood  of  deep  signal  fades.  Inspections  of  beam  images  at  the  receiver-plane  disclose  a  consistent 
pattern;  for  moderate  levels  of  uncorrected  turbulence,  beam  jitter  is  the  dominant  effect. 

The  examples  shown  in  Figure  5-3  provide  some  insight  into  the  mechanisms  that  govern  beam 
statistics  when  the  Wq//q  ratio  is  of  order  unity.  The  far-field  profile  is  generally  characterized  by  a 

single  main  lobe  that  is  somewhat  larger  than  a  diffraction-limited  beam,  which  becomes  more  distorted 
as  vv'o/rQ  increases.  In  this  regime,  however,  irradiance  fluctuations  at  the  receiver’s  boresight  are  driven 

primarily  by  beam  motion  if  tilt  correction  is  not  implemented.  Note  that  the  small  speckle  features  that 
were  prominent  in  Figure  5- 1  are  not  found  in  this  image. 
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RADIAL  DISPLACEMENT  (m) 

Figure  5-3 .  Far-field  beam  image  (left  figure)  and  horizontal-axis  profile  (right  figure)  for  wq//q  =1  . 
The  Hufnagei-  Valley  turbulence  profile  was  used  in  this  simulation  of  vertical  propagation,  which 
resulted  in  an  r{)  value  of  19  cm  at  1.55  mm.  The  propagation  distance  is  40  Mm,  and  the  linear  dimension 
of  the  image  field  is  1  km. 
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Within  the  adaptive-optics  community,  calculations  of  the  far-field  Strehl  are  usually  derived  from 
estimates  of  the  accumulated  phase  errors  over  the  entire  propagation  path.  Since  the  engagement 
geometries  typically  involve  targets  that  are  well  outside  of  the  atmosphere,  this  approach  is  equivalent  to 
collapsing  all  of  the  channel  perturbations  into  a  single  phase  screen  within  the  near  field  of  the 
transmitter.  For  Kolmogorov  turbulence,  Noll  has  shown  that  the  phase  variance  due  to  atmospheric  tilt 
jitter  is  nearly  an  order  of  magnitude  larger  than  the  combined  contributions  of  the  higher-order  Zernike 
modes  [11].  A  scintillation  analysis  based  solely  on  beam  jitter  provides  significant  insight  into  the 
underlying  physical  mechanisms  that  dominate  receiver  fluctuations  when  vv0  «  r0  . 


A  rough  estimate  of  fluctuation  statistics  due  to  turbulence-induced  tilt  can  be  made  by  assuming 
that  the  far-field  beam  profile  and  the  tilt-jitter  distribution  are  both  described  by  Gaussian  functions. 
Using  p  and  0  to  represent  the  radial  and  vectorial  angles  with  respect  to  boresight,  the  normalized  beam 
profile  at  range  L  can  be  written  as 


/(/?)*  exp 


(  n  \ 

2' 

-2 

P 

kwlJl) 

(5-12) 


where  w,  is  the  far-field  beam  radius;  the  distribution  of  the  beam  displacement  angle  would  be 
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(5-13) 


where  at  is  the  single-axis  jitter  radius.  The  tilt  distribution  can  be  converted  to  an  on-axis  Strehl 
distribution  by  combining  the  last  two  expressions 
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where 
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{wJLf 

4  a,2 


(5-15) 


Similar  relationships  have  been  previously  published  by  Fried  and  Titterton  [17,18].  Simulation  results 
for  uncorrected  turbulence  support  this  model,  as  discussed  in  Section  6. 
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An  estimate  of  the  power-law  parameter,  m,  is  obtainable  given  knowledge  of  wt  and  a, .  The  first 
parameter  follows  from  Eq.  (4-9) 


2  L 

w,  = - 

k  w0 


(5-16) 


and  the  second  derives  from  a  calculation  of  the  Zernike  tilt  component  of  turbulence  developed  by 
Sasiela  [8] 
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,  2  5/3  1/3 
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which  yields  the  result 


m  «  0. 196 
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V  r0  J 
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The  scintillation  index  obtains  from  the  moments  of  Eq.  (5-14) 
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so  that 


erf  eJ-4-1  = 


m[m+  2) 


(5-19) 


(5-20) 


1 0/3  5/3 

Note  that  this  predicts  that  m  will  be  proportional  to  (  vvq/jq)  ~  when  the  ratio  is  small  and  ( H'o/ro ) 
when  it  is  large,  which  is  in  qualitative  agreement  with  prior  estimates  [19,20]. 

Unfortunately,  the  analytical  approach  just  described  provides  little  guidance  with  regard  to  the 
scintillation  properties  of  a  tilt-corrected  transmitter.  For  that  case,  one  would  expect  that  irradiance 
fluctuations  will  be  primarily  driven  by  the  changing  dimensions  of  the  beam’s  main  lobe  and  variations 
in  the  combined  strength  of  the  secondary  lobes.  On  the  other  hand,  Eq.  (5-20)  does  suggest  that  the 
scintillation  index  should  be  roughly  proportional  to  the  transmitter-plane  phase  variance,  which 
according  to  Noll’s  estimates  are  all  proportional  to  the  5/3  power  of  the  h'q/tq  ratio  [11] 
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total  phase  variance: 


(5-21) 


<ji*lM(D/rofi=5M(w0/r0f} 


tilt-removed  phase  variance: 


aj*  0.134  (D/r0f13 


0.758  (H.0/r0)5/3 


(5-22) 


Following  this  reasoning,  one  might  expect  the  functional  forms  for  the  uncorrected  and  tilt-corrected 
scintillation  index  to  be  similar,  but  with  leading  constants  that  differ  by  about  an  order  of  magnitude. 


(n  Figure  5-4,  the  scintillation  index  results  from  the  POPS  code  that  were  introduced  previously 

fj  u 

are  compared  with  models  of  the  form  cjj  =  a  (wq/j-q)  within  the  transition  region  between  the  point- 

source  and  large-aperture  limits.  As  before,  the  coherence  diameter  is  computed  for  focused  beam 
propagation  from  0  to  L,  as  defined  in  Eq.  (3-1) 


rn  = 


0.423  k~ 


J(l-z/L)5/3C„2(z)A 


-3/5 


:  focused  beam  . 


In  the  regime  where  near-field  phase  perturbation  dominate  the  shape  of  the  far-field  beam,  the 
scintillation  index  function  displays  a  strong  power-law  behavior.  A  least-square  fit  to  the  uncorrected 
turbulence  data  yields  the  result 

L 

aj  &  1.30  (wq//q)5^  =  0.550  k1  J(l-z/Z,)5/^  C„  (z)  dz  ,  (5-23) 

0 


which  is  consistent  with  predictions  based  on  Eq.  (5-20).  The  shape  of  the  tilt-removed  scintillation  index 
is  similar,  but  a  7/3  power-law  exponent  provides  a  better  fit  to  the  data  than  the  5/3  exponent  used  in  the 
last  expression 
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(5-24) 


These  results  are  judged  to  be  highly  significant,  in  that  the  same  set  of  expressions  have  been 
successfully  applied  to  both  an  uplink  and  a  horizontal-path  geometry  for  which  the  Cn  profile 
descriptions  are  very  different.  Additional  simulation  runs  are  planned  to  assess  the  range  of  applicability 
of  these  two  relationships. 
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Figure  5-4.  Comparison  of  simulation  results  (red  and  blue  data  points)  with  the  transition-region 
scintillation  models  (dashed  lines).  The  upper  figure  represents  an  uplink  geometry  and  the  lower  is 
horizontal  propagation  to  10  km. 
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5.3  SCINTILLATION  INDEX  FOR  w0/r0  »  I :  SATURATION  DUE  TO  STRONG  NEAR¬ 
FIELD  TURBULENCE 


In  the  limit  of  far-field  beam  saturation,  the  irradiance  profile  is  characterized  by  multiple  lobes  of 
roughly  equal  intensity  separated  by  deep  nulls.  The  angular  dimensions  of  both  the  bright  and  dim 
regions  are  of  the  order  of  \jk  .  Two  representative  images  are  shown  in  Figure  5-5. 


RADIAL  DISPLACEMENT  (m) 


Figure  5-5.  Far-field  beam  image  (left  figure)  and  horizontal-axis  profile  (right  figure)  for  wq  /  /  q  =5. 
The  Hiifnagel-Valley  turbulence  profile  was  used  in  this  simulation  of  vertical  propagation ,  which 
resulted  in  an  r()  value  of  19  cm  at  1.55  mm.  The  propagation  distance  is  40  Mm,  and  the  linear  dimension 
of  the  image  fields  is  1  km.  The  average  dimension  of  the  individual  speckles  is  approximately  diffraction- 
limited. 


Superficially,  this  picture  bears  a  resemblance  to  the  beam  profile  for  the  small-beam  case  shown  in 
Figure  5-1,  but  the  physical  processes  involved  are  very  different.  Recall  that  for  the  point-source  limit, 
the  long-exposure  beam  width  is  essentially  diffraction-limited,  and  the  superimposed  modulation  pattern 
is  characterized  by  relatively  shallow  spatial  features  having  an  angular  extent  that  is  small  compared  to 
T/vvq  .  In  the  saturation  limit,  the  size  of  the  long-exposure  profile  is  of  the  order  of  T/aq  ,  and  each  of  the 
speckle  features  is  essentially  an  independent  diffraction-limited  beam  with  random  intensity  and  offset 
angle.  Since  the  superposition  of  the  ensemble  of  beams  is  coherent,  the  null  regions  between  the  peaks 
tend  to  be  very  deep. 
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Although  the  physical  structure  of  this  chaotic  behavior  seems  very  complex,  the  statistical 
descriptions  of  the  fields  in  both  the  transmitter  and  receiver  planes  are  surprisingly  simple.  An  important 
result  of  the  numerical  simulations  performed  as  part  of  this  investigation  is  the  discovery  that  far-field 
irradiance  saturation  is  strongly  correlated  with  saturation  of  phase  errors  in  the  near  field  of  the 
transmitter.  If  the  optical  phase  at  each  point  in  the  transmitter  plane  is  described  by  the  residual,  <9, 
within  the  range  ±n  such  that  ^  =  *9  +  /?  2x ,  saturation  occurs  when  the  distribution  of  9  becomes 
uniform 


p(3)  =  — — : -7i  <  9  <  n  .  (5-25) 

2k 

This  condition  is  clearly  at  odds  with  the  Rytov  approximation,  which  assumes  that  the  magnitude  of  *9  is 
small  compared  to  k and  that  its  distribution  is  Gaussian. 

One  way  to  model  the  far-field  effects  of  the  phase  distribution  given  in  the  previous  expression  is 
to  divide  the  transmitted  beam  into  an  ensemble  of  N  uncorrelated  sources.  It  is  reasonable  to  expect  that 

the  characteristic  spatial  dimension  of  each  source  will  be  of  the  order  of  /q  ,  so  that  N  cc  (  vvq / )  .  If 
each  discrete  source  is  described  by  the  phase,  9( ,  in  the  Fourier  transform  plane,  the  on-axis  Strehl  is  the 
phasor  sum  [21] 
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1  ^ 

2>xp{i>,} 


(5-26) 


If  the  distribution  of  the  random  variable,  9t ,  is  uniform,  the  distribution  of  S'  is  specified  by  the  negative 
exponential  function 

/?(S)  =  N  exp{-N  5}  ,  (5-27) 


which  is  confirmed  by  the  numerical  simulations.  The  first  two  moments  of  this  model  and  the  associated 
scintillation  index  are  summarized  in  the  following  expressions: 
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(5-28) 


While  this  scintillation  index  derived  from  this  model  is  not  sensitive  to  the  magnitude  of  N  (as  long  as  N 
is  large),  an  estimate  of  this  number  is  required  to  compute  the  mean  Strehl  ratio.  The  earlier  prediction 

that  N oc(vv’o/Vq)  is  discussed  in  Section  5.5,  which  provides  an  overview  of  the  Strehl  models 
investigated  in  this  study. 
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5.4  COMBINING  THE  SCINTILLATION  INDEX  MODELS 


Arguments  have  been  presented  in  this  section  to  support  the  concept  that  beamwave  scintillation 
can  be  divided  into  at  least  three  distinct  regions  of  behavior — the  small-aperture  geometry  in  which 
phase  perturbations  are  concentrated  in  the  far-field  of  the  transmitter,  the  wq  *  scenario  for  which 

phase  perturbations  are  concentrated  in  the  near-field  of  the  transmitter,  and  saturated  scintillation  due  to 
strong  perturbations  in  the  near-field  of  the  transmitter.  The  second  case  has  been  subdivided  into 
uncorrected  and  tilt-corrected  (figure)  turbulence;  the  latter  case  being  the  more  desirable  configuration 
and  the  one  more  likely  to  be  implemented  in  a  practical  system.  Since  beam-jitter  correction  is  never 
perfect,  these  two  examples  provide  upper  and  lower  limits  on  system  performance. 

The  scintillation  index  models  that  have  been  developed  thus  far  can  be  summarized  as  follows: 


far-field  effects  that  are  accurately  modeled  by  the  Rytov  approximation . 
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near -fie  Id  effects  driven  by  tilt  and  figure  distortions: 
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and  saturation  due  to  near- field  perturbations: 


(5-32) 


The  following  heuristic  expressions  combine  these  relationships  in  such  a  way  that  a  smooth  transition 
between  each  of  the  regions  is  accomplished 
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Recall  that  Aq  is  defined  for  the  focused-beam  geometry 


/h  = 


1-3/5 


0.423  k2  |(l-z/Z.)5/3  Cl(z)dz 


The  exponential  term  in  the  denominators  of  Eqs.  (5-33)  and  (5-34)  introduces  a  small  amount  of 
overshoot  at  the  initiation  of  saturation,  which  is  consistently  observed  in  the  results  of  the  numerical 
simulations. 

It  has  been  suggested  by  some  authors  that  systems  designed  to  operate  near  the  minimum  of  the 
scintillation  index  curve  provide  optimal  performance.  Since  channel  conditions  can  be  highly  variable, 
this  strategy  may  not  be  practical,  but  it  is  worth  noting  that  this  design  point  can  be  computed  from  the 
composite  <37  expressions  given  above 


uncorrected  turbulence  minimum; 


«2.78£-2 


\z2(\-zll)~'ftCl(z)dz 

0 _ 

J(l-2, !Lfncl(z)dz 

0 


(5-35) 


tilt-corrected  turbulence  minimum: 
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A  comparison  of  the  composite  scintillation  index  models  derived  in  this  section  and  simulation 
results  obtained  from  Lincoln  Laboratory’s  POPS  code  are  shown  in  Figure  5-6.  The  agreement  between 
the  models  and  data  is  very  good,  and  as  noted  before,  the  same  pair  of  expressions  has  been  applied  to 
both  the  uplink  and  horizontal-path  propagation  geometries. 
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Figure  5-6.  Comparison  of  simulation  results  (red  and  blue  data  points)  with  the  composite  scintillation 
index  models  (dashed  lines).  The  upper  figure  represents  an  uplink  geometry  and  the  lower  is  horizontal 
propagation  to  JO  km. 
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5.5  STREHL  RATIO  MODELS 


Investigations  of  scintillation  often  fail  to  include  any  discussion  of  the  average  Strehl  ratio,  which 
is  surprising  in  thst  accurate  knowledge  of  this  quantity  is  essential  for  any  system  application  where  the 
optical  source  is  limited  by  size,  weight,  or  electrical  power  constraints.  This  section  provides  Strehl 
approximations  for  uncorrected  and  tilt-corrected  turbulence;  both  relationships  are  shown  to  be  highly 
reliable  over  a  wide  range  of  channel  conditions. 

The  extended  Marechal  approximation,  which  relates  the  Strehl  to  the  transmitter-plane  phase 
variance,  is  used  extensively  by  the  adaptive-optics  community 

(S)«exp{-o-Jj  .  (5-37) 


9 

To  estimate  the  value  of  the  Zernike-component  estimates  given  in  Eqs.  (5-21)  and  (5-22)  are 
typically  referenced.  While  this  approach  is  usually  adequate  when  the  phase  variance  is  small,  alternative 
representations  have  been  developed  that  are  accurate  for  a  wider  range  of  a ^  values.  A  notable  example 
is  one  published  by  Andrews  and  Phillips  [13] 
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A  number  of  authors,  including  Yura  [7],  have  drawn  a  distinction  between  the  short-exposure  tilt- 
included  and  short-exposure  tilt-removed  turbulence  coherence  diameter.  An  understanding  of  this 
distinction  is  very  important  for  systems  that  incorporate  high-bandwidth  closed-loop  beam  tracking.  In 
the  scintillation  index  plots  shown  earlier,  it  was  demonstrated  that  tilt  effects  become  a  dominant  concern 
when  hjo/Vq  is  of  order  unity.  It  is  not  surprising,  therefore,  that  differences  between  the  uncorrected  and 

tilt-corrected  Strehl  ratios  would  become  significant  in  this  region.  The  following  relationship,  which  is  a 
modified  version  of  Eq.  (5-38),  eliminates  the  tilt  component  of  the  Zemike  phase  variance  until 
saturation  is  achieved 


+ 


5.56- 


4.84 


1  +  0.04  (h'o/'o) 


,5/3 


(w’o/ro) 


<5/3 


-6/5 


(5-39) 


Comparisons  of  these  models  with  results  obtained  from  the  POPS  code  are  shown  in  Figure  5-7.  These 
last  two  relationships,  in  conjunction  with  the  expressions  for  the  scintillation  index  given  previously  in 
Eqs.  (5-33)  and  (5-34),  provide  a  highly  accurate  means  of  estimating  the  low-order  statistical  moments 
of  the  receiver-plane  irradiance  for  a  wide  range  of  propagation  conditions. 


40 


io1  r 


10°  r- 


'  *  *  *  'I 


A 

\  ' 


1  1  lul 

•  Uncorrected  Beam 
a  Tilt-Corrected  Beam 
-  -  -  -  Uncorrected  Model 
- Tilt-Corrected  Model 


£  io-1 

H 

C/5 

<  10’2 

LU 


N  \ 


vertical  geometry 
X  =  1.55  pm 


10'3  r  L  =  40x1 06  m 


r  = 19  cm 
o 


10" 


1  ■  1  1  1 1 J 1  1  1 


10 


10 


\\ 

V 


4, 


X  .  . 


10° 


io’ 


10 


io1  r 


i  i  i  i  i  i  i  1 1  i  i 


10°  r 


LU  ^0 

DC 

I— 

CO 


I  I  I  I  I  I II I  I  I  I  I  II I 

•  Uncorrected  Beam 
A  Tilt-Corrected  Beam 

-  -  -  Uncorrected  Model 

-  -  -  Tilt-Corrected  Model 


10"  r 


10' 


horizontal  geometry 
X  =  1.55  pm 
L  =  10  km 

C  2  =  1.39x10  16  m'273 


r  =  46.1  cm 

o 


% 

\\ 

\  % 

V, 

\\ 

•v 


10" 


10' 


10'1 


10u 


101 


102 


w  /r 

0  0 


Figure  5-7.  Comparison  of  simulation  results  (red  and  blue  data  points)  with  the  Strehl  ratio  models 
(dashed  lines).  The  upper  figure  represents  an  uplink  geometry  and  the  lower  is  horizontal  propagation  to 
IQ  km. 
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6.  RECEIVER-PLANE  IRRADIANCE  DISTRIBUTION  MODELS 


For  many  system  constructs,  the  estimates  of  the  mean  Strehl  and  scintillation  index  developed  in 
the  last  section  are  adequate  to  provide  a  first-order  estimate  system  performance.  However,  some 
applications,  including  free-space  optical  communications,  require  a  more  detailed  knowledge  of  the 
irradiance  statistics  at  the  plane  of  the  receiver.  In  this  section,  distribution  models  are  developed  for  the 
three  sets  of  conditions  discussed  earlier,  i.e.,  wq/Vq  «  1 ,  u'q//q  a  I  ,  and  wq/aq  »  J  .  The  complexity  of 

the  mathematics  parallels  that  of  the  last  section  in  that  the  models  for  the  large  and  small  aperture  limits 
have  relatively  simple  descriptions,  whereas  the  intermediate  case  is  much  more  difficult  to  accurately 
represent. 

A  successful  scintillation  model  must  satisfy  several  constraints,  in  addition  to  the  obvious 
requirement  for  accuracy.  Mathematical  simplicity  is  a  primary  concern,  since  the  distribution  function  is 
likely  to  become  a  building  block  in  a  high-level  model  of  overall  system  performance.  It  is  also 
important  to  be  able  to  associate  the  free  parameters  of  the  model  with  measurable  properties  of  the 
channel,  such  as  the  Cn  profile.  Statistical  models  that  lack  these  properties  are  of  little  use  in 
engineering  design  studies. 

Before  proceeding,  it  should  be  noted  that  much  of  the  literature  relating  to  this  subject  assumes 
that  log-normal  distribution  is  appropriate  for  almost  all  propagation  geometries.  The  work  of  Andrews  et 
a  I .  [22]  is  a  rare  example  of  an  alternative  statistical  description  that  accurately  describes  scintillation  for 
a  wide  range  of  propagation  geometries.  The  three-parameter  formulation  proposed  by  Andrews  is  known 
as  the  gamma-gamma  model,  which  can  be  described  as  the  product  of  a  power-law  function  and  the 
K-Bessel  function 


p{s)= 


2  (a  P) 


[a+0)/2  /  \(a+P)/2-l 


r(a)r(/))(s)((s). 


K 


(6-1) 


The  a  and  /?  parameters  can  be  related  to  the  strength  of  large-  and  small-scale  phase  perturbations, 
respectively.  Although  the  models  developed  in  this  section  have  a  somewhat  different  structure,  they 
share  many  of  the  essential  properties  of  the  gamma-gamma  function. 

6.1  DISTRIBUTION  MODEL  FOR  w0/r0  « 1 :  WEAK  FAR-FIELD  TURBULENCE 

Recall  from  Section  2  that  the  Rytov  approximation  assumes  a  multiplicative  field  perturbation  of 
the  form 


ue  =  exP  {x  +  i</>} 


(6-2) 
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which  indicates  that  irradiance  fluctuations  derive  from  the  real  part  of  the  argument,  Citing  the  central 
limit  theorem,  it  can  be  argued  that  both  x  an^  <t>  are  Gaussian  random  variables,  in  which  case  the 
associated  Strehl  distribution  will  be  log-normal 


p(s)~ 


-exp 

5 


|n(5)-(ln(S)> 


2  of 


(6-3) 


Since  this  approximation  is  generally  only  applicable  when  the  mean  Strehl  is  close  to  unity,  its  shape  is 
specified  primarily  by  the  scintillation  index,  a]  ■  The  relationship  between  <jj  and  the  turbulence  profile 
in  the  small-aperture  limit  was  specified  earlier  in  Eq,  (5-4). 

A  representative  set  of  comparisons  between  simulation  results  and  the  log-normal  model  is  shown 
in  Figure  6-1  for  a  case  where  the  beam  diameter  is  small  compared  to  the  coherence  diameter 
(vvo/rO  -0-0744).  The  shape  of  the  model  functions  is  solely  determined  by  the  value  of  a/,  which  is 
essentially  identical  for  the  uncorrected  and  tilt-corrected  measurements  in  this  example. 

6.2  DISTRIBUTION  MODEL  FOR  vv0/r0  *  1  :  MODERATE  NEAR-FIELD  TURBULENCE 


This  case  is  the  most  difficult  of  the  three  to  model  due  to  the  complex  phenomenology  in  the 
transition  region  between  the  point-source  and  large-beam  limits.  For  uncorrected  turbulence,  tilt-jitter 
becomes  a  dominant  effect,  and  beam-breakup  strongly  influences  the  received  signal  fluctuations  in 
systems  that  incorporate  closed-loop  track  correction.  For  both  of  these  system  architectures,  the  shape  of 
the  irradiance  distribution  is  found  to  be  very  different  from  the  log-normal  function  that  derives  from  the 
weak-turbulence  approximation. 

The  essential  constituent  of  the  tilt-dominated  distribution  is  the  power-law  model,  which  was 
introduced  in  Eq.  (5-14).  Unfortunately,  this  model  does  not  adequately  capture  irradiance  behavior  near 
the  distribution  mean,  because  it  does  not  include  the  effects  of  high  spatial-frequency  phase  fluctuations 
in  the  transmitter  plane.  The  approach  taken  in  this  section  parallels  the  design  of  most  adaptive-optics 
phase  compensation  systems  in  that  the  tilt  and  figure  (tilt-removed  phase)  components  of  the  near-field 
distortion  are  addressed  separately. 
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Figure  6-1.  Comparison  of  uplink  simulation  results  for  Wq  /  /q  = 0.0744  (black  data  points)  with  the  log¬ 
normal  distribution  model  (red  line).  The  upper  figure  represents  the  tilt-corrected  turbulence  case  and 
the  lower  is  the  uncorrected  case. 
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The  analytical  treatment  developed  in  this  section  relies  heavily  on  a  variation  of  the  Nakagami-m 
distribution  that  is  frequently  used  to  characterize  RF  signal  fading  due  to  multipath  [23,24].  The  more 
familiar  form  of  this  expression  describes  the  statistics  of  the  field  amplitude,  but  when  written  in  terms  of 
the  optical  irradiance,  the  following  relationship  is  obtained 


m 


r(m)(s){m(s) 


s 


vW-l 


exp 


(6-4) 


This  representation,  which  is  a  variation  of  the  gamma  distribution  [25],  includes  a  power-law  function 
that  is  dominant  when  S  is  small  and  an  exponential  term  that  defines  the  shape  of  the  curve  near  the 
distribution  mean.  However,  it  is  not  possible  to  separately  adjust  the  properties  of  these  two  functional 
components  because  both  are  constrained  by  the  parameter,  m. 

A  three-parameter  model  version  of  the  Nakagami  model  was  developed  as  part  of  this 
investigation.  The  goal  was  to  retain  the  simplicity  of  the  model,  while  introducing  a  third  parameter  that 
allows  the  shapes  of  the  curves  near  the  mean  and  the  low-Strehl  tail  to  be  independently  controlled.  The 
derivation  begins  with  the  assumption  that  Eq.  (5-14)  can  be  generalized  as  follows 


p{s)  =  < 


f(s/s0) 

0 


:  0  <  S  <  S0 

:  elsewhere 


(6-5) 


where  Sq  <  1  is  the  peak  Strehl  achievable  when  tracking  errors  are  perfectly  corrected.  If  the  tilt- 
removed  Strehl  parameter  is  treated  as  a  time-varying  random  variable,  then  Eq.  (6-5)  can  be  written  as 
the  marginal  distribution 


p(S)=  lp(S;S0)p(S0)dS0  =  mS""'  JW(S0)<i% 


m 


(  V«-l  00 


s 


(So){(S») 


\p{S0)dS0 


J  s 


(6-6) 


The  approximate  form  of  Eq.  (6-6)  is  valid  if  the  standard  deviation,  oq,  of  the  variable  Sq  is  small 
compared  to  its  mean.  This  results  in  a  representation  that  is  the  product  of  distinct  tilt  and  figure 
components. 

In  the  previous  section,  it  was  demonstrated  that  weak  far-field  perturbations  yield  scintillation 
statistics  that  are  accurately  described  by  a  log-normal  distribution.  However,  simulations  show  that  as  the 
perturbations  in  the  near  field  become  stronger  the  statistics  of  the  tilt-removed  Strehl  more  closely 
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resemble  a  simple  Gaussian.  For  a  Gaussian  distribution,  the  integral  indicated  in  the  previous  expression 
is  equivalent  to  the  complementary  error  function 


]P(50)^0=ierfc(^M}  .  (6-7) 

Unfortunately,  the  error  function  is  difficult  to  manipulate  analytically,  so  it  would  be  desirable  to 
substitute  a  more  tractable  expression  having  comparable  properties  near  the  distribution  mean.  From  a 
Taylor’s  series  expansion  centered  on  (Sq),  it  can  be  shown  that  the  complementary  error  function 

closely  resembles  the  exponential  expression 


—  erfc 
2 


\Sq-(So)\ 

1  JlcQ  J 


!n(2) 


IW 


y' 

) 


where 


J2 Tn_{Sp) 

!n(2)  £70 


(6-8) 


This  is  equivalent  to  replacing  the  Gaussian  distribution  with  a  form  of  the  Weibull  distribution  [25], 
specifically 


yf2ft 


exp 


crn 


(Sq-(Sq)) 

2al 


Sa 


ln(2  )fi 

(So)  {(So) 


exp 


-ln(2) 


Sn 


\0 


(So). 


(6-9) 


The  observations  outlined  above  motivated  the  formulation  of  the  following  three-parameter 
scintillation  model,  which  is  referred  to  in  this  report  as  the  modified  Nakagami  (or  mod-N)  distribution 


p(S)  = 


C/3 


,  m-\ 


r(mip) (s){  {s) 


c- 


exp 


-  C 


(S) 


(6-10) 


r(lm  +  l)/0) 

The  constant  C  =  — - — - - - — -  has  been  defined  to  simplify  the  notation.  Note  that  this  result  is 

equivalent  to  Eq.  (6-4)  when  /3  =  1,  but  that  non-unity  values  of  /?  provide  a  means  to  independently 
control  the  characteristics  of  the  tilt-removed  component  of  the  expression.  The  moments  of  this  function 
have  a  simple  analytical  form 


Mk 


r((m+k)//i) 

fMT 

r{m/p) 

v c ) 

(6-11) 


so  the  scintillation  index  can  be  expressed  as 
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(6-12) 


j2  r({m+ 2)/p)r(m/fi)  t 

r((m  +  l)/j3f 

The  integral  of  the  probability  density  is  the  incomplete  gamma  function 


S 


P{S)=  \P{S')dS' 


0 


r[mlP,(CSl(S)f} 


(6-13) 


Although  the  proposed  model  and  its  associated  moments  have  a  relatively  simple  mathematical  structure, 
a  strategy  for  computing  the  model  parameters  m  and  ft  is  not  immediately  obvious.  The  following  steps 
outline  a  straightforward  process  for  developing  estimates  of  both  parameters  given  knowledge  of  the 
scintillation  indices  for  uncorrected  and  tilt-corrected  turbulence. 


If  one  accepts  the  premise  that  the  effects  of  the  figure  component  of  near-field  phase  errors  can  be 
modeled  as  a  Weibull  distribution,  then  the  following  result  is  derived 


Pfeis o)  =  - 


r(i//i)(r(i//3)  s0 


(So)  [  P  (So) 


\0-[ 

exp  < 

) 

~ 

nm  h 

P  (S0). 


\0 


(6-14) 


Note  that  this  equation  does  not  include  the  power-law  parameter,  m.  Therefore,  the  scintillation  index  for 
the  figure  component  is 


r(\+i!P) 

r(m//?)2 


(6-15) 


where  the  approximation  given  in  the  second  relationship  is  accurate  when  ft  >  1  .  This  step  establishes 
the  value  of  ft.  An  estimate  of  the  power-law  parameter,  m,  derives  from  Eq.  (6-12) 


2  r((m  +  2)/p)r(m/fi) 
r({m  +  \)/p)2 


^r(1+01^) 

Jj 


(6-16) 


where  the  approximate  form  provides  an  acceptable  estimate  when  mj ft  <  1 .  If  required,  the  accuracy  of 
the  estimates  for  the  two  shape  parameters  can  be  improved  by  inverting  Eqs.  (6-15)  and  (6-16) 
numerically. 


The  following  three  figures  illustrate  the  evolution  of  the  tilt-corrected  and  uncorrected  irradiance 
distribution  functions  for  increasing  values  of  Wq/tq.  All  simulations  incorporate  the  Hufnagel-Valley 

i 

C„  profile  and  vertical  propagation.  The  log-normal  model  is  included  for  comparison. 
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Figure  6-2.  Comparison  of  uplink  simulation  results  for  h'q/zq  =  0.187  ( black  data  points)  with  the  log¬ 
normal  distribution  model  (red  line).  The  upper  figure  represents  the  tilt-corrected  turbulence  case  and 
the  lower  is  the  uncorrected  case. 
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Figure  6-3.  Comparison  of  uplink  simulation  results  for  wq/tq  =  0.371  (black  data  points)  with  the  log¬ 
normal  distribution  model  (red  line).  The  upper  figure  represents  the  tilt-corrected  turbulence  case  and 
the  lower  is  the  uncorrected  case. 
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STREHL 

Figure  6-4.  Comparison  of  uplink  simulation  results  for  wq  /  /  q  =0.930  (black  data  points)  with  the  log¬ 
normal  distribution  model  (red  line).  The  upper  figure  represents  the  tilt-corrected  turbulence  case  and 
the  lower  is  the  uncorrected  case. 
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The  previous  examples  demonstrate  that  the  modified  Nakagami  model  provides  a  very  good  fit  to 
the  uncorrected  turbulence  cases,  but  is  less  successful  in  matching  the  shapes  of  the  distributions  for  tilt- 
corrected  turbulence  when  vvq / / q  is  small.  In  the  latter  case,  the  log-normal  distribution  is  clearly 

superior,  presumably  because  the  relative  contribution  of  far-field  turbulence  is  significant.  A  method  for 
combining  the  log-normal  and  modified  Nakagmi  functions  is  described  in  Section  6.4. 

6.3  DISTRIBUTION  MODEL  FOR  w0/r0  »  1 :  SATURATION  DUE  TO  STRONG  NEAR¬ 
FIELD  TURBULENCE 

In  the  discussions  relating  to  saturated  scintillation  presented  in  Section  5,  it  was  noted  that  this 
condition  appears  to  be  associated  with  near-field  distortions  that  are  strong  enough  to  create  a  uniform 
distribution  of  phase  errors  within  the  range  of  ±k  .  In  this  limit,  the  field  in  the  transmitter  plane  can  be 
modeled  as  an  ensemble  of  N  independent  emitters,  where  N  ~1 j(S).  The  resulting  distribution  is  the 
single-parameter  negative-exponential  function 

5 

(S> 

This  is  a  special  case  of  the  Weibull  distribution  given  in  Eq.  (6-14)  for  which  the  shape  parameter,  /?,  is 
equal  to  1,  and  is  also  a  special  case  of  the  modified  Nakagami  function  with  m  /?  1  .  A  distinguishing 

characteristic  of  negative  exponential  distribution  is  that  the  function  approaches  a  constant  for  small 
values  of  S. 

The  scintillation  index  for  this  distribution  is  unity,  and  the  shape  of  the  curve  is  entirely 
determined  by  the  mean  Strehl,  (5).  In  the  saturation  limit,  the  mean  Strehl  can  be  estimated  from  Eq. 
(5-38) 

(5)  «  0.125  (w0/r0)  2  ,  (6-18) 

which  indicates  that  the  number  of  equivalent  incoherent  emitters  is  proportional  to  the  number  of 
coherence-diameter  regions  within  the  profile  of  the  transmitted  beam. 

A  representative  example  of  the  irradiance  distribution  for  saturated  scintillation  is  given  in  Figure 
6-5.  The  shapes  of  the  curves  for  tilt-corrected  and  uncorrected  turbulence  are  nearly  identical  and  are 
solely  determined  by  the  value  of  the  mean  Strehl.  It  is  interesting  to  note,  however,  that  the  uncorrected- 
turbulence  plots  in  both  Figure  6-4  and  6-5  show  evidence  of  a  slightly-negative  power  law  exponent  (i.e., 
m  <  l ).  This  condition  appears  to  coincide  with  the  region  of  overshoot  in  the  scintillation  index  function 
shown  in  Figure  5-6  and  may  be  the  result  of  beam  tilt  fluctuations  that  exceed  the  short-exposure  beam 
radius. 


Af)=^exp< 
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Figure  6-5.  Comparison  of  uplink  simulation  results  for  r§-  3.72  (black  data  points)  with  (he 
negative-exponential  distribution  model  (red  line).  The  upper  figure  represents  the  tilt-corrected 
turbulence  case  and  the  lower  is  the  uncorrected  case. 
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6.4  COMBINING  THE  IRRADIANCE  DISTRIBUTION  MODELS 


The  statistical  models  discussed  in  the  last  two  sections  are  appropriate  for  disturbances  confined 
primarily  to  the  near-field  of  the  transmitter,  and  they  all  represent  variations  of  the  three-parameter 
function  that  has  been  referred  to  in  this  report  as  the  modified  Nakagami  distribution.  Far- field  phase 
perturbations,  on  the  other  hand,  result  in  signal  distributions  that  closely  resemble  a  two-parameter  log¬ 
normal  function.  A  formalism  that  achieves  a  smooth  transition  between  these  two  dissimilar 
representations  is  presented  in  this  section. 

A  review  of  the  tilt-corrected  simulation  data  that  appears  in  Figures  6-2  and  6-3  provide  a  good 
starting  point  for  this  discussion.  In  the  first  example,  the  data  are  closely  matched  to  the  log-normal 
function,  whereas  in  the  second  case  the  simulation  results  in  the  tail  of  the  distribution  fall  between  the 
estimates  provided  by  the  near-  and  far-field  scintillation  models.  In  both  cases,  the  modified  Nakagami 
distribution  overestimates  the  likelihood  of  a  deep  fade,  although  it  could  be  argued  that  this  discrepancy 
is  typically  not  important  when  the  scintillation  index  is  small.  However,  for  communication  systems 
designed  to  operate  near  the  link  margin,  a  higher  degree  of  prediction  accuracy  may  be  needed. 


The  approach  proposed  for  combining  the  near-  and  far-field  scintillation  models  derives  from  the 
discussion  of  the  composite  scintillation  index  provided  in  Section  5.5.  The  following  expression 
represents  a  weighted  average  of  the  two  distributions 


p(s)= 


J  FF 


Pff  {S)  + 
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J  NF 


.  Pnf(s) 
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FF 


C J i 


-]2 


J  NF 


(6-19) 


where  the  values  of  the  weights  are  specified  by  the  scintillation  indices  as  defined  in  Eqs.  (5-29)  through 
(5-32).  The  far-  and  near-field  functions,  ppp  and  pNF  ,  refer  to  the  log-normal  and  modified  Nakagami 
distributions  defined  in  Eqs.  (6-3)  and  (6-10).  The  linear  addition  indicated  in  the  last  expression 
preserves  both  the  statistical  mean  and  variance  of  the  two  component  distributions. 

The  effect  of  the  composite  distribution  model  is  most  significant  in  the  region  where  the  far-field 
and  near-field  contributions  are  comparable.  In  Figure  6-6,  new  plots  of  the  tilt-removed  data  shown 
earlier  in  Figures  6-2  and  6-3  are  given  along  with  the  composite  model  shown  in  the  dashed  green  curve. 
In  the  former  case,  the  composite  model  essentially  replicates  the  log-normal  model  since  far-field  effects 
dominate,  but  in  the  second  example  it  provides  a  significantly  better  match  to  the  simulation  output  than 
either  of  the  component  distribution  functions. 


54 


> 

H 

c n 


LLl 

o 

> 


DO 

< 

CO 

O 

DC 

a. 


£ 

O) 

z 

Ui 

D 

> 


CD 

< 

CD 

O 

DC 

a 


STREHL 


Figure  6-6.  Comparison  of  uplink  simulation  results  for  the  tilt-corrected  beam  ( black  data  points)  with 
the  composite  distribution  model  (green  line).  The  upper  figure  represents  the  wq/tq  -0.187  case  and 

the  lower  is  the  result  for  wq  //q  -  0.371 . 
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7.  CONCLUSIONS 


Numerical  simulations  of  beamwave  propagation  frequently  reveal  significant  discrepancies 
between  the  scintillation  properties  of  the  computational  output  and  predictions  based  on  the  Rytov 
theory.  This  problem  has  frequently  been  encountered  in  studies  of  engagement  geometries  of  practical 
interest,  including  many  that  are  relevant  to  free-space  optical  communications.  System  designers  are  now 
requesting  models  that  more  accurately  predict  channel  behavior,  but  which  are  also  simple  enough  to 
allow  their  integration  into  high-level  engineering  software.  The  relationships  developed  in  this  report 
represent  an  attempt  to  satisfy  both  of  those  requirements. 

Analytical  investigations  of  weak  scintillation  have  generally  been  restricted  to  two  cases  that  can 
be  solved  through  the  use  of  the  Rytov  approximation — propagation  of  an  infinite  plane  wave  and 
illumination  from  a  point  source.  Neither  of  these  geometries  is  an  adequate  representation  of  a 
transceiver  link  that  transmits  a  beam  of  finite  dimensions,  and  which  imparts  a  focus  that  is  matched  to 
the  receiver  range.  Both  of  these  propagation  constraints  are  fundamental  to  most  systems  of  practical 
interest,  and  both  have  been  applied  in  this  analysis. 

One  of  the  primary  goals  of  this  study  is  a  much  better  understanding  of  the  physical  processes  that 
drive  receiver-plane  signal  fluctuations.  Two  parameters  are  found  to  be  of  particular  importance.  The 
first  is  the  ratio  between  the  outgoing  beam  diameter  and  the  turbulence  coherence  diameter;  this  result  is 
not  unexpected  in  that  this  ratio  appears  frequently  in  publications  dealing  with  turbulence  effects  and 
adaptive  optics.  The  second  factor,  which  is  somewhat  more  subtle,  is  the  relationship  between  locations 
of  the  individual  perturbation  planes  and  the  Rayleigh  range  of  the  transmitted  beam.  This  study  has 
identified  a  fundamental  difference  between  the  effects  of  phase  distortions  introduced  in  the  near  field  of 
the  transmitter  (within  the  Rayleigh  range)  and  the  far  field  (beyond  the  Rayleigh  range).  In  general,  near- 
field  perturbations  are  not  well  characterized  by  predictions  based  on  the  Rytov  approximation.  This 
observation  is  particularly  important  for  applications  such  as  uplink  communications,  which  exhibit  a 
strong  near-field  component. 

This  report  provides  a  practical  approach  to  predicting  irradiance  statistics  given  knowledge  of  the 
beam  size,  propagation  geometry,  and  turbulence  profile.  Although  the  relationships  derived  have  thus  far 
only  been  validated  for  the  cases  of  vertical  uplink  and  horizontal  propagation,  it  is  believed  that  they  can 
be  accurately  applied  to  a  much  wider  range  of  engagement  geometries. 
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APPENDIX 


The  following  derivation  represents  an  alternative  approach  to  the  estimation  of  the  scintillation 
index  due  to  far- field  perturbations  and  is  intended  to  support  the  result  developed  in  Section  5.1.  As  in 
the  previous  discussion,  it  is  assumed  that  the  primary  effect  of  phase  distortions  close  to  the  receiver 
plane  is  to  modify  the  local  curvature  of  the  propagating  wavefront.  In  essence,  the  physical  model  is  that 
the  beam  propagates  through  a  succession  of  concave  and  convex  lenses  before  being  intercepted  by  the 
receiver. 

To  obtain  a  quantitative  description  of  the  irradiance  changes  that  result  from  the  focus-driven 
model  proposed  above,  the  focus  perturbation  term,  £*,  will  first  be  defined  such  that  the  effective  field 
curvature  in  the  z/  plane  is 

1  1  /,  X 

-  =  -(,+£)  .  <A-I) 

J I  Z1 


The  perturbed  beam  diameter  in  the  Z2  plane  further  along  the  propagation  path  is 


d  7  —  d 9  +  d 


~7,  “  \{z2  ~z\)  =  d2+£~(z2-z\)  , 

\f\  Z\J  z\ 


(A-2) 


and  the  modified  boresight  irradiance  is 

In 


r2  = 


1  +  £  d\  (z2  —  Z|  )/c/2  z\  ^Jr£{z2~z\)lz7 


/2  [l-2£- (z2 -Z|  )/z2  . 


(A-3) 


At  this  point  the  assumption  will  be  made  that  s  is  a  zero-mean  Gaussian  random  variable  with 

9 

variance  .  If  the  value  of  this  parameter  can  be  related  to  the  turbulence  strength  at  z/,  then  the 

statistics  of  I2  are  easily  derived  from  Equation  (A-3).  Noll’s  article  [ll]  shows  that  within  a  circular 
aperture  of  diameter  d ,  the  focus  component  of  the  phase  error  is 


focus 


=  0.023 


VOJ 


=  0.0097  d5p  k 2  C2n  (z)  Az 


(A-4) 


where  Az  is  the  phase-screen  depth.  For  a  wavefront  characterized  by  the  curvature  radius,/  the  piston 
and  tilt-removed  phase  variance  within  the  circle  is  [8] 
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( A-5) 


<y  a 
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In  this  context  the  parameter,  y  -  1  //  ,  can  be  regarded  as  a  statistical  quantity  with  expectation  value 


<.->= 


768 

77 
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If  x  =  1  f  f  is  a  zero-mean  Gaussian  variable,  then  the  distribution  for_y  is  the  gamma  function 


p{y) 


exp 


which  has  the  mean  value,  ^  y)  =  c2x  .  This  provides  the  key  relationship  between  the  variance  of  1/  / 

and  the  focus  component  of  the  phase  variance  given  in  Equation  (A-6).  When  this  result  is  applied  to  the 
perturbation  term  in  Equation  (A-l)  and  combined  with  the  expression  for  the  characteristic  radius 
developed  in  Section  4.3 
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it  can  be  seen  that 
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The  final  step  of  the  derivation  incorporates  the  relationship  between  I2  and  8  given  in  Equation  (A-3)  to 
develop  an  expression  for  the  receiver-plane  intensity  variance 
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}  *  2.64  A' 7/6  z?/6 


Z2 


-z. 
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( A-9) 


The  variance  for  the  entire  propagation  path  derives  from  the  assumption  that  the  phase  perturbations  are 
additive  and  independent 
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C„2(z)dz 


(A- 10) 


a}  =  2.64  it7/6  p6 


P5/6 


v  z2  ; 


The  leading  constant  in  this  result  differs  from  the  standard  expression  given  in  Equation  (5-3)  by  about 
20%,  but  this  qualitative  agreement  tends  to  support  the  argument  that  the  primary  effect  of  far-field 
perturbations  is  embodied  in  a  model  that  describes  far-field  perturbations  as  a  series  of  thin  lenses. 
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